arXiv:1504.05525v4 [nlin.CD] 14 Dec 2016 


Delocalization of Disturbances and the Stability of AC Electricity Grids 


Stefan Kettemann 

Jacobs University, Department of Physics and Earth Science, Campus Ring 1, 28759 Bremen, Germany and 
Division of Advanced Materials Science Pohang University of Science and Technology (POSTECH) San 31, 

Hyoja-dong, Nam-gu, Pohang 790-78Jf., South Korea 
(Dated: December 15, 2016) 

In order to study how local disturbances affect the AC grid stability, we start from nonlinear 
power balance equations and map them to complex linear wave equations. Having obtained sta¬ 
tionary solutions with phases (pi at generator and consumer nodes i, we next study the dynamics of 
deviations. Starting with an initially localized perturbation, we find it to spread in a periodic grid 
diffusively throughout the grid and give the parametric dependence of diffusion constant D. We 
apply the same solution strategy to general grid topologies and analyse their stability against local 
perturbations. The perturbation remains either localized or becomes delocalized, depending on grid 
topology, power capacity and distribution of consumers and generator power Pi. Delocalization is 
found to increase the lifetime of perturbations and thereby their influence on grid stability, while 
localization results in an exponentially fast decay of perturbations at all grid sites. These results 
may therefore lead to new strategies to control the stability of electricity grids. 


The stability of electricity grids requires to protect 
them against perturbations [THS]. Therefore, electrical 
power systems must be constructed in such a way that a 
physical disturbance does not result in exceeding bounds 
of system variable fluctuations. The energy transition 
towards an increased supply of decentralized renewable 
energy necessitates to study consequences of such struc¬ 
tural changes for the stability of electricity grids and to 
find efficient ways to modify them to ensure their sta¬ 
bility [4]. Since this is a highly complex and nonlinear 
problem the study of its dependence on network topol¬ 
ogy, operating conditions and forms of disturbances re¬ 
quires to make modeling assumptions [2]. Recently, the 
synchronization of rotor angles in electricity grids has 
been modeled by a network of nonlinear oscillators [5]-[7]. 
Here, networks of generators and engines are described 
by a system of coupled differential equations for local ro¬ 
tor angles of generators and loads pi, where i are grid 
nodes. The numerical solution of these differential equa¬ 
tions showed that, on the one hand, networks become 
more unstable with increasing decentralization against 
perturbations on short time scales with large amplitude, 
while the danger of a blackout can be reduced by de¬ 
centralization [6]. In this article, we study how phase 
perturbations evolve with time in AC grids. The origin 
of such phase perturbations may arise for example due to 
local fluctuations in generating power from a wind gen¬ 
erator. We start by finding stationary solutions for the 
spatial distribution of phase pi, for given distribution of 
active and reactive power. Pi and Qi at the grid nodes 
i. Next, we reconsider the nonlinear dynamic power bal¬ 
ance equations. For small deviations from the stationary 
solutions, we derive linear wave equations describing the 
phase perturbation dynamics. Solving these equations, 
we explore how a local phase perturbation propagates 
with time through the grid. Depending on the geograph¬ 
ical distribution of power, grid power capacity and topol¬ 
ogy we find that it may either spread diffusively or be¬ 
come localized. This phenomenon is generally known as 


Anderson localization [8], where the coherent scattering 
of waves in a random medium causes their localization. 

Steady state power flow in AC transmission grids .— 
The power balance equations in AC transmission grids 
are obtained from Kirchhoff’s laws at node i as 


Si = Y,Vi 
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where Si = Pi P \Qi. Pi is the active power produced 
at generator nodes Pi > 0, or consumed at consumer 
nodes < 0, satisfying total power balance condition 
= 0 with N the total number of nodes. Vi is 
the voltage at node i. In an AC grid the reactive power 
Qi of consumers is given while the one of generators is 
adjusted[3l |9]. The transmission line from node i to j 
has impedance Zij. Neglecting small losses due to Joule’s 
heat, we assume them to be purely inductive Zij = luoLij, 
where Lij is the transmission line inductance between 
nodes i and j, uj the grid frequency. Then, the voltage 
at node i is Vi = V exp(i(/?^), where V is fixed to nominal 
grid voltage and the power capacity of a transmission 
line is Kij = V‘^ j{ujLij)Aij, where Aij is the adjacency 
matrix of the grid. Note that only N — 1 phase angles pi 
remain to be determined as function of the distribution 
of power Pi at N nodes constrained by the total power 
balance condition, so that reactive power Qi is fixed at 
all nodes. Defining = ex.p{—ip^{t)), with Pi{t) = 

ujt P 0^, phase angle 0^ at node i in steady state, it is 
convenient to write Eq. 0 as a linear wave equation, 

sp (t) = Yj Pi (*) - (i)) • (2) 
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Phase dynamies .— Phase dynamics in AC electricity 
grids has been modeled by active power balance equa¬ 
tions with additional terms describing the dynamics of 
rotating machines [5HZl EHIl]- One term describes the 
inertia to changes of kinetic energy of a synchronous ro¬ 
tating generator or motor with rotor angle pi with inertia 
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J, when we assume that all loads are either synchronous 
or induction motors, whose dynamics can be modeled 
that way [3]. Another term describes the damping with 
coefficient 7 . Adding these terms to the active power bal¬ 
ance equations, the real part of Eq. ( [T|), yields for purely 
inductive transmission lines, I3E1I31T1HI1] 

We note that for fixed voltage V there are no dynamic 
terms in the reactive power balance equation, which ap¬ 
pear only in higher order when also voltage dynamics in 
addition to the phase dynamics is considered [15]. 

Dynamics of Disturbances in the Grid. — In order to 
study the propagation of disturbances, such as fluctua¬ 
tions in power supply dP(t) or in power capacitance 6Kij^ 
we set (fi (t) = ujt ai (t) with steady state phases 

the solutions of Eq. ([^. We study the dynamics of 
phase disturbances ai{t) which are governed by 

dfUi + 2Tdtai = ^ - E ^ sin(0° - 9° + Ui - Q:j),(4) 



FIG. 1: Left: Square grid of generators (red) and con¬ 

sumers (blue). Arrows: direction of active power trans¬ 
mission F. Right: Topology of german transmission 
(380kV(red), 220kV(blue)) and high voltage distribution grid 
(110kV(black))[T6]. A propagating disturbance is sketched as 
red circles. 


where k is a wave number. Eor periodic boundary condi¬ 
tions on the grid of linear size L in all d directions with 
unit vectors e^, n = l,...,d, '^/.(x + Le^) = ' 0 /c(x), for 
all n = 1 ,..., d, the wavenumber is k = 27TnlL^ where the 
components of the vector n, Un^n = l,...,d are integers. 
Eor each wave number k we find a solution of the form 
Eq. where the phase factors fjQk i 'f’Ck are given by 


where F = ^jj. Considering small perturbations from 
the stationary state we expand Eq. Q in 
yielding linear wave equations on the grid, 

+ 2Tdtai = ~ ( 5 ) 

j 

with hopping amplitude tij = Kij cos{0^ — 0^)/{Juj). 
Note that tij depends both on power capacitance Kij 
and thereby the grid topology, as well as on the initial 
distribution of power Pi through the stationary phases 
6 >?. With ai{t) = X]„c„iexp(-ia;„t) we get for T = 0 
the spectral representation of the linear wave equation 

^n^ni ~ ^ ^ ^ij {^ni ^nj\ (b) 

j 

with eigenfrequencies uj^ and eigenmodes c^. cq = 0, 
cJo = 0 correspond to the stationary solution. Eor F 7 ^ 
0 we get from Eq. (§ the same eigenmodes with 
complex eigenfrequencies fin = — lE + i^/F^ — ujI Eq. 

can be solved for arbitrary electricity grids like the 
german transmission grid shown in Eig. 1 (right), where 
the hopping matrix elements tij with stationary phases 
are obtained from the solution of Eq. Q. 

Square grid. — As an example, let us first consider a 
grid where all transmission lines have equal length a. We 
start with a periodic arrangement of generators and con¬ 
sumers as shown for a square grid in Eig. 0(1 eft). Assum¬ 
ing that all generators on sublattice G generate power 
Px = X G G, while all consumers on sublattice C 
consume power Px = — P, x G G, we find the solution by 
making the Bloch Ansatz for the stationary solution 

'ipk{^eC,t) = 'ipcke'^'^exp{-iLot), (7) 


i’Ck = exp (i4) \l)Gk with sin4 = P/(fkK), ( 8 ) 

where fk = 2 Yln=i cos(/c^a) depends on wave number 
components kn- Eor given reactive power Q, which is 
for this arrangement of consumers and generators con¬ 
strained by Eq. to be the same at all nodes, the wave 
vector k is determined by the equation 

fl = {Q/K-2df+PyK\ (9) 


The transmitted power between sites i and j is Fij = 
iKij{fji{t) — 'ipj{t))'ipi{ty. Eor the homogenous state 
k = 0 the active power transmitted between neighbored 
sites i G G and sites j G G is ReFij = P/(2d), as shown 
in Eig. [^(left) on a square lattice, where arrows indicate 
the direction of active power transmitted from genera¬ 
tors to neighbored consumers. The reactive power in the 

transmission lines is ImF = K{1 — — 4 ^?^), which 

is for P ^ P much smaller than the transmitted ac¬ 
tive power. Eor finite wave vector k the active power 
transmitted between i G G and j G G for k = kCn is 
ReFij = K sm{6k ± knu), with unit vector between i and 
j, Cij = Pen. In all other directions ReFij = K sinSk. 

Stability .— Disregarding the dependence of ai on the 
perturbation at neighbored sites aj reduces Eq. Q to 
the one of a damped, driven nonlinear pendulum. Eor 
large times t ^ 0 it has two stable solutions: 1. Sta¬ 
tionary solution dtOLi = 0 , = n27r, n integer, to which 

small deviations decay exponentially. 2. Overswinging 
pendulum solution, when driving force and damping are 
in balance. The phase velocity converges then to 


j 


Kij sin(f2if + 9° - 9^ + rji) 

jpWTIP 


( 10 ) 









3 



FIG. 2: Phase curves of Eq. 0 for OLi^ setting otj — 0. 
Flow directions are indicated by arrows. Stable fixed points 
are oli — n27r, dtoa — 0. Separatrices are an = 2Sk + (2/ + 
l)7r, dtai = 0. Blue shaded areas are basins of attraction. All 
other phase points converge to Eq. (10) (red curve). 


with frequency fti = Pi/{2rjuj) and phase shift r]i = 
arctan(2r/^2i). In Fig. 2 phase curves are shown for 
the square lattice with 0^ — 0^ = + knCi. Regions 

of stability are shaded in blue, outside which all phase 
points converge to open orbit solution Eq. (10) (red). 
A condition for phase points to lie inside the stability 
region is obtained by approximating its irregular shape 
by an ellipse whose vertex on the a-axis is given by the 
saddle point, an = 26^ — tt. Its vertex on the dta-ax.is 
is given by a linear interpolation of stable trajectory on 
the separatrix dta = —(tt — 2 J/c)tanC, where tan^ = 

are inside the basin of attraction if they satisfy 


* ^ tan2 Ci ^ ■ 


( 11 ) 


Propagation of Local Disturbances .— Now, we can 
study the propagation of disturbances in a square grid 
by inserting the analytical steady state solutions Eqs. 
dM ) into the dynamical equations Eq. 0 . If we per¬ 
turb the phase at node i = an{t = 0) = ao 0 , while 
^^( 0 ) = 0 at all other nodes i ^ n., that disturbance 
excites nodes i 7 ^ n at later times t > 0. For a small 
perturbation, satisfying Eq. ( pAj ), we can use linearized 
wave equations Eq. ^ with tij = Kij cos 6k/(Jou). Us¬ 
ing the spectral representation ai{t) = ^ 

insertion into Eq. 0 gives the complex frequency 


eq = -iF 




2dK cos 5k 
Jwr2 


(1 
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For finite F and small q the dispersion is quadratic, 
~ — q^. For large momenta the ballistic 

limit, the dispersion is linear eg|r=o = with velocity 



FIG. 3: (color online) Phase curves at initial site ri and other 
sites rn , . Time progresses as indicated in color bar in units 

of T = 1/F. Delocalization causes slow decay of the perturba¬ 
tion at the initial site, and an initial increase, followed by a 
slow decay at other sites (upper curves). Localization causes 
an exponentially fast decay at all sites (lower curve). 


Vk = y /K COS 6k / {Juj)a. An initially localized perturba¬ 
tion, Cq = const, becomes for times t > r = 1/F and dis¬ 
tances exceeding the mean free path I = vr^ |g — n| > ^ 


ai{t) = 


0^0 


{fLKDktj aYl^ 


exp 


ADkt 


(13) 


Thus, the initially localized perturbation spreads diffu¬ 
sively with diffusion constant 


Dk = 


Ka^ cos 6k 

UJJ 


(14) 


as shown in Fig. 3, where dai is plotted versus phase 
perturbation ai at initial site i = I and other sites i'^i" 
for diffusive times t > r = 1/F, using Eq. (13). Time 
is progressing as indicated on the color bar in units of 
r. Diffusion causes slow decay of the perturbation at 
the initial site, and an initial increase, followed by a slow 
decay at other sites. The area At where nodes become 
perturbed at time t, increases diffusively with time as 
At = Dkt. In order to quantify the stability with Eq. 

Jill we calculate the spatial average, 6a = ^^ a^ = 

decaying 

with a power law in time. For power capacity K = .5G1U, 
d = 2, grid frequency uj = 2ii 50/s, power P = 1.9G1U, 
inertia J = 10^ kgm^ and damping F = 1/s the perturba¬ 
tion spreads initially with velocity v = 2.22ajs. Beyond 
mean free path I = 2.22a it spreads diffusively with dif¬ 
fusion constant = 2.46a^/s. The time to diffuse a 
length L = 10a is tn = L‘^/Dq = 40.64s. The closer 
power P is to capacity limit Kc = 2dK, the smaller is 
diffusion constant the longer survives the disturbance. 

The propagation of any type of disturbance can be 
studied with Eq. 0. F. e. a change in power 6P at to = 
0 at neighbored nodes i,j, 6Pi = —6Pj = 6P changes 
transmitted power between nodes /c, / at time t. 


6Fki{t) = ±6 PA 


kl 


TT^a^ 

uioDkC 
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with + when /c is a generator and I a consumer and — 
vice versa. Similarly, the change in power flow due to a 
static perturbation SPi = —SPj = SP is given by 


SFki = ±SPAki2TTa'^/{ri - rj)^. 


(16) 


An instantaneous change in power capacitance 5KAij 
causes a change of power Eq. ( p!^ , and Eq. (16) for 
a static perturbation, replacing 6P by 6K, respectively. 

Localization of Disturbances. — Eq. can be applied 
to study the propagation of disturbances in AC grids of 
arbitrary distribution of power P^, arbitrary grid topol¬ 
ogy and power capacity Kij by calculating the hopping 
amplitude tij = Kij cos{0^ — 0^)/{Juj) from the steady 
state solutions as found by solving Eq. In a 

first attempt to model the complexity of a real grid, we 
take a random distribution of tij caused by the wide dis¬ 
tribution of power Pi and the real grid topology with 
its complex network structure. Eq. first appeared 
in the problem of randomly coupled atoms in harmonic 
approximation. Eor chains it has been solved for var¬ 
ious random distributions of tij nHHii]- If tij is taken 
from a box distribution, the density of eigenmodes is 
constant, p(fln) = 1/^? for 0 < For nonzero 

On the eigenstates are localized with localization length 
^ 1/On [IHH 2 T]. The localization length of the 
lowest eigenfrequency Oi ^ uja/L is large, of the order 
of the system size ^ L. The highest eigenfrequency 
On ^ cj has a localization length ^ ^ 2a, twice the 
length a of a transmission line. In dimension d = 2, which 
is the situation most relevant for real grids, all eigen¬ 
states remain localized for nonzero On, but the localiza¬ 
tion length can be for small On exponentially large, so 
that it typically is larger than the system size for realistic 
grid extensions L. This might explain that in Ref. m we 
found in square grids a long range power law decay with 
power g = 2 as in Eq. (16), even when taking a random 


distribution of Pi. Also in a real grid topology we found 
long range decay [17]. Typically, the localization length 
is smallest in treelike grids, while it becomes larger, the 
more meshed the grid becomes. In dimensions d > 2, 
there is a critical value cjc, such that for ujn > ujc ah 
modes are localized, while they are extended for uJn < uJc 
[ 2 TJ [ 22 ] . If the phase perturbation is initially in a state 


localized around site ro with localization length it 

decays in time t as aAt) = exp(— \ exp(— 

^ _ sn 

where = Re\r — — cj^]. The typical phase and 

frequency shift is found to decay exponentially fast as 


5a = ao\/^n/<^exp(-fnt), 5uj = Vn5a. (17) 


Thus, localization causes exponentially fast decay of 
phase perturbations at all nodes i, as shown in Eig. 3, 
where it is compared with a delocalized phase perturba¬ 
tion decaying with a power in time t, Eq. (13). 


Conclusions .— Local perturbations, arising for exam¬ 
ple from power fluctuations, are found to spread diffu¬ 
sively in a periodic grid, decaying slowly with a power 
law in time and space. The closer the generator power P 
comes to the capacity limit Kc^ the smaller is diffusion 
constant D and the longer takes the perturbation to de¬ 
cay. Modeling the complexity of a realistic grid with a 
random distribution of generators and consumers and/or 
random transmission power capacity the phase perturba¬ 
tion is found to become localized in 1- and 2-dimensional 
grids. Localization leads to an exponentially fast decay 
of phase perturbations at all sites, while delocalization 
results in diffusive slow decay, Eig. Initially small 
perturbations may then add up at some nodes to large 
perturbations and push the system outside of the region 
of stability. We conclude that it is favorable for sta¬ 
ble grid operation to ensure that phase perturbations re¬ 
main localized, decaying exponentially fast at all sites. 
The consequences of these results for real el ectr icity grid 
topologies will be studied by solving Eqs. numeri¬ 

cally in a future publication. We also plan to study how 
the spreading of perturbations is modified when includ¬ 
ing voltage fluctuations, using for example the third order 
model [Il[3l[7|. While we assumed here a network of syn¬ 
chronous generators and motors, modern wind turbines 
are rather induction generators, the most modern ones 
being the doubly fed induction generator, converting the 
power from AC to DC and then to AC with the grid 
frequency [ 3 ]. Thus, the energy transition towards an in¬ 
creased supply of decentralized renewable energy neces¬ 
sitates to get a better understanding how the dynamic 
equations are modified and to understand the resulting 
consequences for the grid dynamics. 
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